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Abstract 

This study presents two new modeling strategies. First, a methodology for 
representing the physical process of stormflow within a TOPMODEL framework is 
developed. In using this approach, discharge at quickflow time scales is simulated and a 
fuller depiction of hydrologic activity is brought about. Discharge of water from the 
vadose zone is permitted in a physically realistic manner without a priori assumption of 
the level within the soil column at which stormflow saturation can take place. 
Determination of the stormflow contribution to discharge is made using the equation for 
groundwater flow. No new parameters are needed. Instead, regions of near saturation 
that develop during storm events, producing vertical recharge, are allowed to contribute 
to soil column discharge. These stormflow contributions to river runoff, as for 
groundwater flow contributions, are a function of catchment topography and local 
hydraulic conductivity at the depth of these regions of near saturation. The second 
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approach improves groundwater flow response through a reduction of porosity and field 
capacity with depth in the soil column. Large storm events are better captured and a 
more dynamic water table develops with application of this modified soil column profile 
(MSCP). The MSCP predominantly reflects soil depth differences in upland and lowland 
regions of a watershed. Combined, these two approaches — stormflow and the MSCP — 
provide a more accurate representation of the time scales at which soil column discharge 
responds and a more complete depiction of hydrologic activity. Storm events large and 
small are better simulated, and some of the biases previously evident in TOPMODEL 
simulations are reduced. 

GAP Index Term Numbers: 1800, 1836, 1860 

Introduction 

Successful modeling of the hydrologic cycle requires representation and 
quantification of the various pathways by which water migrates through a catchment. 
Over the last twenty years, the capabilities and fidelity of hydrology models have 
increased considerably. Many factors impacting the spatial and temporal variability of 
the hydrologic cycle, and the scales at which these processes operate, have been 
elucidated and applied to model simulations. These developments have permitted 
modeling of the spatial distribution of soil moisture levels across the land surface, water 
movement within the soil column, and the relative contributions of evapotranspiration 
and river runoff to the efflux of water from the catchment system. 




Many hydrology models have been developed using TOPMODEL formulations 
(Ambroise et al., 1996a; Beven and Kirkby, 1979; Koster et al., 2000; Sivapalan et al., 
1987; Stieglitz et al., 1997). A framework of analytic equations, TOPMODEL is based 
on the idea that topography is the primary determinant of the distribution of soil moisture 
at and within the land surface (Beven, 1986a; Beven, 1986b; Beven and Kirkby, 1979; 
Beven et al., 1994). 

TOPMODEL formulations define areas of hydrological similarity — that is, points 
within a watershed that respond to meteorological forcing in similar fashion, saturating to 
the same extent, producing the same levels of discharge, etc. These points of 
hydrological similarity are identified by an index that is derived from analysis of 
catchment topography. This topographic index is often of the form ln(a/tanfi), where 
tanfi is the local slope angle at a patch on the land surface, and a is the amount of upslope 
area draining through that patch. Lowland areas tend toward higher topographic index 
values, due to a combination of either low slope angle or large upslope area. Upland 
areas tend conversely toward lower topographic index values. Points within a catchment 
with the same topographic index value are assumed to respond identically to atmospheric 
forcing. Thus within a TOPMODEL framework, the topographic index provides the 
fundamental unit of hydrological response. 

This fundamental topographic unit is derived from three basic assumptions (see 
Ambroise et al.. 1996a; Beven, 1997 for details):(l) the water table is approximately 
parallel to the topographic surface so that the local hydraulic gradient is close to tanfi\ (2) 
the saturated hydraulic conductivity falls off exponentially with depth; and (3) the water 
table is recharged at a spatially uniform, steady rate that is slow enough, relative to the 
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response time scale of the watershed, to allow the assumption of a water table distribution 
that is always at equilibrium. These assumptions permit reconstruction of the spatial 
variability of catchment response to meteorological forcing solely from modeling of the 
response of the mean state. This quasi-stochastic approach is at once computationally 
efficient while still permitting dynamic representations of physical processes within the 
system. 

These conceptual underpinnings provide a foundation for physically-based 
modeling of catchment hydrology. Land surface models based on these formulations 
have been applied to catchments both large (Duchame et al., 2000) and small (Beven and 
Kirkby, 1979). However, while the TOPMODEL framework has provided hydrologists 
with a powerful and efficient tool for modeling hydrologic conditions, the full 
complement of dynamic processes is not represented in most TOPMODEL applications. 
As a consequence, model simulations often perform poorly during drier conditions. The 
response of watersheds to wetting by spring snowmelt and to storms after an extended 
dry period have proven particularly difficult to represent. Correction of these simulation 
inaccuracies requires a more detailed depiction of the hydrologic cycle and catchment 
physical structure. 

Discharge response 

Four physical processes contribute to river runoff in a watershed: 1 ) precipitation 
onto stream channels; 2) overland flow; 3) shallow subsurface storm flow; and 4) 
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groundwater flow (Homberger et al., 1 998). The first two of these processes respond 
very rapidly, producing spikes in hydrographs during and immediately after storm events. 

The third mechanism, shallow subsurface stormflow, responds at the quickflow 
time scale. (N.b. In this paper the terms stormflow and groundwater flow refer to 
physical processes; the terms quickflow and baseflow define short and long time scales of 
response.) Tracer studies have shown that shallow’, subsurface regions of the soil column 
can support significant levels of flow during storm events (Dewalle and Pionke, 1994, 
Hendershot et al., 1992; Ogunkoya and Jenkins, 1993). Such regions can exist as perched 
water tables, disconnected from the true w ; ater table supporting groundwater flow 7 (Gile, 
1958; Hammermeister et al, 1982; Noguchi et al, 1999; Wilson et al., 1990). These 
perched water tables, by virtue of their development in the vadose zone nearer the land 
surface, can flow more quickly, discharging their waters more rapidly to the catchment 
river network. However, the timing of the development of these perched water tables and 
their size and placement within the soil column have proven difficult to model. 

Stormflow has been represented in Variable Infiltration Capacity (VIC) models 
(Liang et al., 1994; Lohmann et al., 1998a; Lohmann et al., 1998b). VIC models provide 
a viable modeling alternative to the more physically-based TOPMODEL approach. 

These parameterized reservoir models mimic the time scales of catchment hydrologic 
response to storm events and can be structured to emulate quickflow and baseflow 7 time 
scales. Such time scales are calibrated to inferred rates of quickflow 7 and baseflow 
derived from runoff analyses. 

Recently, Scanlon et al. (2000) included stormflow in a TOPMODEL-based 
hydrology model through introduction of a parameterized quickflow reservoir. This 



model represents a hybridization of the parameterized reservoir and TOPMODEL 
approaches. The authors used a saturation deficit model as their basic soil column 
structure, then added a second saturation deficit reservoir from which stormflow would 
be determined. The groundwater flow and stormflow components were partitioned so as 
not to overlap and thus supersaturate the soil column. Analyses of hydrograph response 
to different storm events and antecedent conditions, and piezometer sampling of the study 
catchment were used to determine the saturation deficit recession coefficients. Two time 
scales of response, an order of magnitude apart — 471 and 36 hours for groundwater flow 
and stormflow respectively — were delineated; these rates were applied to the reservoir 
discharge formulations. The authors used a linear rate of decrease in the transmissivity of 
the soil column for the stormflow reservoir, while maintaining an exponential decay of 
transmissivity for the groundwater reservoir. Flow from the subsurface stormflow 
reservoir to the groundwater reservoir was facilitated by a linear recharge function. For 
the study watershed, the model depicted both the rapid and slow discharge of water from 
the soil column following storm events. 

This hybrid modeling approach (Scanlon et al., 2000) produced three additional 
parameters: one for stormflow reservoir maximum capacity; a second for the stormflow’ 
discharge rate; and a third for the linear recharge function. Successful simulation of 
catchment hydrology required calibration of these parameters to rates inferred from 
analysis at a highly instrumented experimental watershed. An attractive alternative to 
this hybrid model approach w'ould be one allowing a more general inclusion of 
stormflow, w’hile requiring fewer parameterizations. Indeed, if stormflow could be 
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incorporated in a manner consistent with the TOPMODEL formulations that govern the 
flow of water within the soil column, no additional parameters would be necessary. 

The fourth process contributing to river runoff in a watershed, groundwater flow, 
provides most of the baseflow during extended periods between storm events. This 
discharge mechanism is represented in land surface hydrology models constructed in the 
TOPMODEL framework (Ambroise et al., 1996a; Ambroise et al., 1996b; Stieglitz et al., 
1997; Wood et al., 1990). Specifically, TOPMODEL formulations permit dynamically 
consistent calculations of both the partial contributing area — from which precipitation 
onto stream channels and overland flow can be determined — and the groundwater flow 
that supports this area. However, model recession hydrographs often run high during wet 
periods, then low during drier months, as will be shown in this paper. Consequently, 
calibrations applied during wet conditions corrupt simulation accuracy during dry 
conditions, and vice versa. 

These model discharge biases reflect inappropriate levels of groundwater flow 
generation. Within the TOPMODEL framework, groundwater discharge is determined 
by the height of the water table (or analogously by the saturation deficit). Low 
groundwater flow simulation during dry conditions stems from greater depth of the water 
table and concomitant lower hydraulic conductivities controlling both discharge and 
recharge rates. During such dry conditions the model water table height often is not 
sufficiently responsive to recharge, and thus groundwater flow discharge is muted. A 
more responsive water table needs to be represented if groundwater discharge rates are to 
be more accurately modeled. 
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Such a dynamic water table could be effected by increasing infiltration-recharge 
rates; however, these rates are set by Darcy's Law and the TOPMODEL assumption that 
saturated hydraulic conductivity decays exponentially with depth — the same attributes 
controlling groundwater discharge. Thus, recharge cannot be altered without changing 
groundwater discharge, or decoupling the mechanisms controlling these two processes 
an unrealistic anisotropy. Instead, a physically representative modification is needed that 
engenders greater response of the water table to recharge, thereby intensifying 
groundwater flow response. 

Here we present two new strategies intended to address these model 
shortcomings: 1) a physically-based approach to modeling stormflow and its application 
in a manner consistent with TOPMODEL assumptions; 2) a modified soil column 
framework, in which porosity and field capacity are realistically allowed to change with 
depth, that provides a more responsive water table and better groundwater flow 
simulation. Together, these strategies produce a more complete and accurate depiction of 
hydrologic activity at the catchment scale. 

The Hydrology Model 

The hydrology model employed for this study has been previously described 
(Stieglitz et al., 1997). Two methods are used for modeling the flow of water within a 
catchment. The first is a soil column model that simulates the vertical movement of 
water and heat within the soil and between the soil surface plus vegetation to the 
atmosphere. The ground scheme consists often soil layers. Layer thicknesses arc 
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structured in a geometric series determined from the depth of the first ground layer 
typically 4 centimeters for this study. Diffusion and a modified tipping bucket model 
govern heat and water flow, respectively. The prognostic variables, heat and water 
content, are updated at each time step. In turn, the fraction of ice and temperature of a 
layer may be determined from these variables. Transpiration and other surface energy 
balance calculations use a standard vegetation model (Pitman et al., 1991) that includes 
bare soil evaporation and canopy interception loss. 

The second method partitions the catchment surface into two distinct hydrologic 
zones: saturated lowlands; and unsaturated uplands. Using the statistics of the 
topography, the horizontal movement of groundwater is tracked from the uplands to the 
lowlands (a TOPMODEL approach). Combining these two approaches (Figure 1) 
produces a three-dimensional picture of soil moisture distribution within a catchment. 
The partitioning of runoff and surface water and energy fluxes is effected without the 
need to explicitly model the landscape. Specifically, an analytic relation, derived from 
TOPMODEL assumptions, exists between the mean water table depth (z), determined 
from the soil column model, and local water table depth at any location x ( z x ) (Sivapalan 
et al., 1987; Wood et al., 1990) 


1 /f\ln(a/tcm{i) 


0 ) 


where ln(a/lan(i) A is the local topographic index at location x; A is the mean watershed 
value of ln(a/tan[i), and / is the rate of decline of the saturated hydraulic conductivity 
with depth in the soil column. By setting z., equal to zero, i.e., locating the local water 
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table depth at the surface, saturated regions of the land surface can be identified. This 
partial contributing area includes all locations for which 


bi( a/tan f3) x > A + yj (2) 


From this partitioning, the contributions to river runoff of both precipitation 
directly onto stream channels and overland flow (saturation-excess runoff) can be 
quantified. Following Sivapalan et al., (1987) groundwater flow (Qb) is: 


a = 


AKfz = 0) 
f 




(3) 


where A is the area of the watershed, and K s is the saturated hydraulic conductivity at the 
surface. This flow through the soil matrix supports river discharge between storm events. 

This combined approach to modeling the land surface has been validated at 
several watersheds, ranging in scale from 2.2 km (Stieglitz et al., 1999) to 570,000 km 
(Duchame et al., 2000). 

The Approach 

A Conceptualization of Siormflow 


Within our hydrology model, Darcian flow accounts for recharge of the water 
table and occurs within layers of the soil column for which volumetric soil moisture is 
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greater than 70% of field capacity (see Hillel, 1 977). Our conceptualization of stormflow 
makes use of this condition. We argue that not only will Darcian flow produce vertical 
recharge of the water table, but it will also bring about stormflow discharge. Such 
stormflow is a natural consequence of the topographic variability of the watershed and 
the TOPMODEL assumption that water tables form parallel to the land surface. While 
our hydrology model conceptualizes a vertical soil column and horizontal land surface, in 
fact most of the land surface is sloped (Figure 2). Thus, gravity will not force Darcian 
flow within the soil column solely in the direction of the water table — a component of the 
flow will be directed laterally. Just as the true water table migrates through the soil 
column, regions of near saturation in the vadose zone will also produce discharge. It is 
thus possible to develop a representation of shallow subsurface stormflow making use of 
the existing formulations for groundwater flow (i.e. Equation 3). As with groundwater 
flow, gravity guides the motions of these stormflow waters. 

Model stormflow is initiated within a layer of the vadose zone when the 
volumetric soil moisture level exceeds 70% of the field capacity (Figure 3). Calculation 
of the stormflow component begins w'ith redistribution of the excess water into a fully 
saturated region from the base of the layer upward: 
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where .vvv, is the stormflow water table in layer /, zb, is the bottom boundary of model 
layer /, 6, is the volumetric soil moisture in layer /, 0?.?, is the field capacity in layer /, and 
0, is the porosity in layer i. As in Stieglitz et al. (1997), groundwater flow is calculated 


given a water table of that height (Q. m , )• Groundwater flow is then calculated a second 
time for a w'ater table of height to the top of the layer below' {Qsi>i)- This second value is 
subtracted from the first, leaving only the flow w'ithin the layer in question. 

Qsi = Qswi - Qsbi ( 5 ) 

Total stormflow for a given time step is merely the sum of the stormflows for each layer 
in the vadose zone. 


Qs = IQsi ( 6 ) 

This formulation requires no new parameterizations. It incorporates the local 
hydraulic conductivities within each layer. Saturated zones that develop nearer to the 
surface flow more quickly than zones formed at depth. This circumstance necessarily 
generates different rates of stormflow and groundwater flow. These differences reflect 
the time scales of quickflow and baseflow. 

A Modified Soil Column Profile (MSCP) 

Catchment upland areas tend to possess shallower soils than lowland areas (Cox 
and McFarlane, 1995; Webb and Burgham, 1997; Yanagisawa and Fujita. 1999). 
Consequently, at the depth of upland bedrock, lowland areas may still hold sedimented 
sands, clays, and organics, supporting a greater porosity and field capacity. This has 
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some important implications. When the catchment is very wet, a given amount of water 
table recharge is spread across the entire catchment, since the water table exists 
everywhere. When the catchment is dry, on the other hand, the water table is limited to 
the lowlands, and a given amount of recharge is allocated to only a fraction of the 
catchment. As a result, in the dry catchment, the level of the water table increases more 
quickly for a given amount of recharge — the dry catchment has a more responsive water 
table. 

A one-dimensional soil column framework, of course, can never allow the explicit 
treatment of this three-dimensional behavior. Nevertheless, a suitable parameterization 
could allow the one-dimensional model to capture this behavior in a gross sense. That is, 
a parameterization could be devised that makes the water table respond more dynamically 
to recharge as the level of the water table decreases. 

We introduce such a parameterization here. Like many hydrology models, our 
starting model used a single parameter value for porosity and a second for field capacity, 
regardless of the depth of the soil column (Stieglitz et al., 1997). With the new 
parameterization, the imposed porosity and field capacity in the soil profile decrease with 
depth into the soil, so that a unit of recharge increases the water table level more quickly 
when the water table is low (that is, when the catchment is dry). Such a modified soil 
column profile (MSCP) is consistent with TOPMODEL steady-state assumptions. This 
MSCP could, in a sense, be interpreted as an "averaging" of soil depth properties over the 
upland and low land areas of the catchment — as the depth into the soil increases, the 
fraction of the catchment with bedrock at that depth increases, so that the average 
porosity decreases. 
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Porosity, in addition, has been shown to decrease with depth at a point, as 
demonstrated in field site analyses using a variety of methods, including mercury 
intrusion (Ajmone-Marsan et al., 1994), bulk density (Asare et al., 1999; Bonell et al., 
1981; Cox and McFarlane, 1995), and fractal approaches (Bartoli et al., 1993; Oleschko 
et al., 2000). These profiles reflect near-surface bioturbation and soil compaction. 
Although the impact of compaction and bioturbation on water table dynamics is probably 
significantly less important than that of spatial variability in bedrock depth, the observed 
one-dimensional porosity profiles are at least consistent with the proposed 
parameterization. 

The MSCP should produce a more responsive water table and generate greater 
groundwater flow. To test this hypothesis, we applied simple and equal geometric 
reductions to both the porosity (cm 3 /cm 3 ) and the field capacity (cm 3 / cm 3 ) - ranging from 
1% to 15% per layer — to the soil column of our hydrology model. By the above 
arguments, steeper catchments with shallower upslope soils should require a greater rate 
of decrease of porosity and field capacity with depth. 

Model Simulations 

Experimental watersheds at Sleepers River and at Black Rock Forest were used in 
this study. These sites are topographically representative of the rolling hillslope and 
steeper ravine catchments that dominate the hydrology of the northeastern United States. 
Meteorological conditions at both watersheds typify the daily, seasonal and interannual 
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variability of this mid-latitude region. Periods of drought, particularly during summer, 
are not uncommon for either watershed. 

At each study site, four model simulations were performed: 1) a baseline run 
without any of the new modifications; 2) simulation with stormflow; 3) simulations with 
the MSCPs; 4) simulations with both stormflow and the MSCPs. 

Results 

Application to the Sleepers River Watershed 

The Sleepers River watershed (111 km 2 ) located in the glaciated highlands of 
Vermont is hydrologically representative of most upland regions in the northeast United 
States. As such, this site was chosen in 1957 as an experimental watershed by the 
Agricultural Research Service (ARS) to provide a better understanding of natural 
watershed behaviour and aid in the development of testing physically based hydrologic 
models (Anderson, 1977). Nested entirely within the Sleepers River Watershed is the W- 
3 sub watershed (8.4 km 2 ). The topography is characterized by rolling hills, and the soils 
are predominantly silty loams. Vegetation cover is approximately equally distributed 
amongst grasses, coniferous forest and deciduous forest. Hourly measurements of air 
temperature, dewpoint temperature, incoming shortwave and thermal radiation and wind 
speed are recorded within the watershed. Average annual air temperature is 4.1°C with a 
standard deviation of 1 1 ,4°C. Mean hourly precipitation is determined from 7 gauges 
placed within the W3-sub-catchment; annual precipitation totals approximately 109cm. 
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Another data set contains the snow water equivalent (SWE), snow depth, snow 
temperature and soil temperature. Total snow depth averages about 254cm and snow 
cover usually persists from December to March with snowmelt in late March and April. 
Hourly runoff data are available from a gauge at which the Pope Brook leaves the W3- 
sub-catchment of the Sleepers River (Stieglitz et al. 1997). Five years of meteorological 
and hydrologic data collected between 1969-74 w'ere used to drive the model. 

Hydrographs of the Sleepers River model runs are shown for years 1971 and 1973 
(Figures 4 and 5). The baseline model runs depict the gross characteristics of runoff 
generation; however, many simulation inaccuracies are apparent (Figures 4a and 5a). 
Specifically, we identify three model shortcomings (see Figure 4a): 

1) Model response to storm events during dry conditions — i.e. October through 
April — is often reduced and spiky, lacking a short-term discharge recession 
curve in the days immediately following a storm. 

2) Groundwater flow response is high during transitions from wet to dry 
conditions — such as the long May through July 1971 recession at Sleepers 
River. 

3) Model response to initial wetting from spring snowmelt is muted and delayed. 

A fourth shortcoming, low bias between storm events during dry periods, is also 
apparent. However, most likely this bias reflects the absence of any depiction of deep 
aquifer water discharge. Such deep flows operate at very long time scales, which are not 
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represented in TOPMODEL dynamics. It is not our intention in this study to address this 
shortcoming. 

With stormflow activated there is improvement in the match of modeled and 
measured hydrographs (Figures 4b and 5b). This model simulation generates river 
discharge at shorter time scales, better representing storm event response, especially in 
dry conditions. Smaller storm events are now resolved, and recession curves after 
summer storm events are evident. There is also considerable improvement in the model 
hydrograph response to spring snowmelt. 

Model simulations with the different MSCPs produce modest improvement in the 
match of modeled and measured hydrographs. A reduction in porosity and field capacity 
of 5% was found to best improve model discharge (Figures 4c and 5c); however, this 
simulation still fails to capture many smaller storm events and for those it does, produces 
little short-term recession in the days immediately following the storm. At the baseflow 
time scale, however, groundwater flow is more responsive, producing greater total 
discharge than the baseline run. Consequently, the May 1971 response is exaggerated, 
but the long July recession following this event is better matched (Figure 4c). 

Sensitivity analyses were performed to determine which MSCP worked best in 
conjunction with stormflow. With both stormflow' and a 5% reduction of porosity and 
field capacity, modeled and measured hydrographs match best (Figures 4d and 5d). Once 
again, the stormflow component generates discharge at quickflow time scales. The 
model is responsive in both wet and dry conditions, simulating discharge recession for 
storms large and small, and best captures catchment response to spring snowmelt. Dry 
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season storm events are well depicted; however, groundwater discharge bias between 
these storms persists. 

Application to the Black Rock Watershed 

The Black Rock forest is a 1500 hectare preserve located in the Hudson 
Highlands region of New York. Elevations in the forest range from 1 10 to 450 m above 
mean sea level, with seasonal temperatures ranging from -2.7 °C to 23.4 °C. The medium 
texture soils are typically very thin, with parent material located from 0.25 to greater than 
lm below the surface in the depressional areas. Soils in the lowland areas are more 
organic than upslope, but bulk densities are not significantly different. Exposed bedrock 
is common throughout the preserve and consequently the area was not extensively farmed 
during the period of European settlement. Lumber extraction ceased in 1927 and the 
forest has been managed as a preserve without significant disturbance since that time. 

The system is typical of the Quercus dominated, secondary growth forests that have 
characterized the NE United States over the past century. The catchment (1.35 km*) is 
drained by a single stream, Cascade Brook. Average hourly discharge from Cascade 
Brook is monitored continuously using a V-notch weir installed in 1998. Hourly 
measurements of precipitation, air temperature, dewpoint temperature, incoming 
shortwave radiation and wind speed are also taken. Hourly thermal radiation for the site 
is calculated following the methodology of Anderson and Baker (1967). Three years of 
meteorological and hydrologic data, collected between 1998-2000, were used to drive the 
model. (Note that the Black Rock weir malfunctioned during November 1999.) 
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Hydrographs of the Black Rock model runs are shown for years 1999 and 2000 
(Figures 6 and 7). The baseline simulations at Black Rock (Figures 6a and 7a) are 
considerably less accurate than their Sleepers River counterparts. This circumstance is 
due in part to differences in catchment structure. The Black Rock watershed is steeper, 
rockier and possesses shallower soils than Sleepers River. Two of the three shortcomings 
apparent in the Sleepers River baseline run are evident from the Black Rock baseline 
simulation: 

1) Dry condition storm events are not depicted. 

2) Groundwater discharge is high during transitions from wet to dry conditions. 

Warm winter conditions precluded significant snowpack development at Black Rock 
from 1998-2000. Consequently, the third shortcoming, poor response to initial wetting 
from spring snowmelt, is absent. Of particular note is the long summer drought of 1999 
that was broken in September by Hurricane Floyd (Figure 6a). The baseline model run 
fails to simulate this flood event, and the smaller storms subsequent. 

The addition of stormflow generates a better discharge simulation (Figures 6b and 
7b). Smaller storm events — particularly during drier periods — are better captured, and 
there is a visible short-term recession curve following the initial flood event. However, 
not all of the characteristics of the measured hydrographs are well matched. In particular, 
large storm events are still not well simulated, and local river runoff maxima, such as the 
overtopping of the weir during Hurricane Floyd in 1999 are still often underrepresented. 
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Model runs with the MSCPs produce considerable improvement of discharge 
simulation. A reduction in porosity and field capacity of 12% was found to best improve 
model discharge (Figures 6c and 7c). The shallow, rocky soils of hillslope and upland 
regions of Black Rock justify this larger rate of decrease in porosity and field capacity. 
The soil profile modification intensifies groundwater flow response to events poorly 
captured in the baseline simulation, producing greater total discharge. The initial flood 
due to Hurricane Floyd is now well depicted, as is the large storm in March 1999 (Figure 
6c). Many smaller storm events, however, remain unresolved, and the MSCP model run 
simulates little of the short-term recession present in the measured hydrographs. 

Analyses were again performed to determine which MSCP worked best in 
conjunction with stormflow. With both stormflow and a 12% reduction of porosity and 
field capacity, modeled and measured hydrographs match best at Black Rock (Figures 6d 
and 7d). Smaller storm events are resolved, and once again the stormflow component 
generates discharge at quickflow time scales. The model is responsive in both wet and 
dry conditions and best captures the intensity of catchment response to large storm 
events. Hurricane Floyd is well captured, as is the large storm in March 1999, and the 
recessions following these flood events are represented (Figure 6d). Some problems 
remain: groundw'ater flow response remains low between storm events during dry 
periods; and runoff generation is too heightened in late 2000. 

Analysis of Water Table Depth 


20 



Further illustration of the impacts of stormflow and the MSCP can be seen from 
examination of modeled, mean water table depth, (D ■ Figure 8 presents comparisons of 
the three new model simulations with the baseline run for the Sleepers River catchment. 

With stormflow activated (Figure 8a) less water recharges the water table than for 
the baseline run. As a consequence, the water table is lower, and groundwater flow is 
necessarily diminished; however, rather than a reduction in hydrograph response, the 
missing water mass is instead shunted to stormflow discharge, raising river runoff levels 
following storm events (Figures 4b and 5b). Not only is the timing of discharge thus 
redistributed, but overall the catchment is more responsive with stormflow activated, 
generating greater total discharge for the five-year duration of the simulation (data not 
shown). 

Model simulation with the MSCP — 5% reduction in porosity and field capacity — 
produces a more dynamic water table (Figure 8b). The modified soil column holds less 
water at depth and therefore is more responsive to an equal volume of recharge. This 
greater responsiveness necessarily intensifies groundwater flow generation and shortens 
the time scales at which this discharge occurs. Additionally, the water table drops more 
abruptly in response to an equal volume of discharge. Thus, the water table rises higher 
in the soil column with wet conditions but lower during dry periods. These changes in 
water table response and groundwater flow production are not offsetting; rather, the 
raising of the water table during wet conditions has a disproportionately larger impact on 
discharge quantities than the lowering during dry periods. Overall, more total discharge 
is generated. This more dynamic water table explains the May 1971 response at Sleeper 


River, in which peak runoff is overly simulated, but July recession is well matched 
(Figure 4c). 

The effects of stormflow and the MSCP on water table depth offset one another 
during wet conditions (Figure 8c). This balance reduces the exaggerated discharge 
simulation found in the model run with the modified soil column profile alone; however, 
the improvement of the July 1971 recession at Sleepers River remains (Figure 4d). 

During dry conditions, the impacts of the two modifications are additive, as both bring 
about a lowering of the water table. As a consequence, underrepresentation of discharge 
persists between dry period storm events. In spite of this bias, water is still available for 
stormflow generation, and the recharge that does reach the water table can raise it more 
quickly, producing faster wetting of the catchment during spring snowmelt. 

Impact of Soil Layer Structure 

Additional analyses were performed to determine the sensitivity of the combined 
stormflow/MSCP simulations to the layer resolution of the model soil column. Layer 
thicknesses were altered by changing the depth of the first ground layer, from which all 
other layer depths are calculated (see Methodology). Runs with both stormflow and the 
MSCP were made at both Sleepers River and Black Rock with first ground layer depths 
of 2. 3, 4. 5. 7 and 12 centimeters. Discharge simulations were found to be insensitive to 
these alterations and therefore appear robust. 


Spectral Analysis 


Whereas the model hydrograph improvements effected by stormflow and the 
MSCP are readily apparent at Black Rock, the impacts at Sleepers River are more 
difficult to discern. To elucidate better the differences at Sleepers River, the faster time 
scales of discharge response were assessed by spectral analysis. Fourier analyses were 
performed on all five years of the measured runoff record and model simulations 
(January 1, 1970 through September 30, 1974). Figure 9 shows the full spectral profile of 
the measured runoff. Most of the power is in the low frequencies, matching the cycle of 
annual snowmelt at Sleepers River. 

Figure 10 presents a partial spectral profile of the higher harmonics (200 and 
higher). The measured runoff record (Figure 10a) has considerably more power than the 
baseline simulation (Figure 10b) at these frequencies. Figure 10c shows the partial 
spectral profile for model simulation with stormflow; power increases over the baseline 
run for harmonics 200 and higher. These frequencies represent responses on time scales 
shorter than nine days, including quickflow time scales. For model simulation with the 
MSCP (Figure lOd), power increases preferentially for harmonics below 400 — time 
scales of four days or longer. Simulation with both stormflow and the MSCP (Figure 
lOe), appears to have an additive effect, combining the power increases of simulation 
with only stormflow' and only the MSCP. This pow'er spectrum best matches the 
magnitudes of the measured runoff record at harmonics above 200, including the highest 
frequencies representing quickflow time scales. These latter time scales reflect discharge 
from intial spring snowmelt and stormflow response. 


Discussion 


Evaluation of Stormflow 

The implementation of stormflow alone produces considerable improvement in 
model representation of catchment discharge. Shallow subsurface stormflow responds to 
storm events on a shorter time scale than baseflow, and thus allows for a more rapid 
flushing of waters from the catchment. This stormflow process is critical during dry 
months. Recharge waters that drain vertically to the deep water table are subject to a 
very low hydraulic conductivity and thus very low rates of groundwater flow. 
Consequently, in the absence of stormflow simulation, storm waters are captured by the 
deep water table and in essence sequestered until wetter conditions prevail, the water 
table rises, and groundwater flow rates increase. This mechanism explains why models 
without stormflow work best in wetter conditions when groundwater flow is more 
responsive. With stormflow activation, however, the catchment is responsive to storm 
events in both wet and dry conditions. This greater responsiveness results from a partial 
redirection of vertical recharge waters — instead of first recharging the water table then 
slowly discharging to river runoff, a portion of the shallow subsurface water moves 
quickly and directly to streambeds. This effect produces both a redistribution of the 
timing of discharge and an increase in total runoff generation. 

While our conception of stormflow is physically based, the perched water tables, 
or fully saturated regions, which can develop above the water table and generate 
stormflow, are not explicitly represented by our modeling approach. Within the land 


24 


surface, these perched water tables can form in several ways. Heterogeneities in the soil 
column — different soil types, even sheets of bedrock — produce surfaces of lower 
permeability upon which water may pool within the soil column. These waters then flow 
downhill over these lower permeable surfaces. This lateral flow continues only as far as 
the edge of the substrate; once this edge is reached vertical recharge of the true water 
table can resume (Noguchi et al., 1999). Because TOPMODEL formulations work from 
the mean, the simulation of such discontinuities in lateral flow is problematic; 
TOPMODEL does not allow for pockets of perched bedrock. Without some type of 
parameterization there is no means of allowing for spatial heterogeneity, i.e. the 
identification of these perched water table surfaces and their edges. 

Another contribution of stormflow to catchment discharge is derived from 
perched water tables that develop in the soil column by virtue of the decay of hydraulic 
conductivity with depth. This decay produces a convergence of water with downward 
flow in the soil column. As waters converge, zones of saturation form. Such zones can 
develop over entire regions of the watershed, allowing for continuous lateral flow to 
discharge areas. Our methodology makes no effort to identify such perched water tables 
explicitly; however, it works in a conceptually similar fashion, using regions of near 
saturation in the vadose zone to generate stormflow. This approach to modeling 
stormflow is also consistent with experimental evidence showing that macropores, 
mesopores, soil pipes, and perched rock within the soil column can direct flow both 
vertically and laterally (Noguchi et al., 1 999; Sidle et al., 2000). While the TOPMODEL 
approach does not detail such flow pathways, our stormflow conceptualization provides 
an implicit representation of such lateral flow in the vadose zone. 


The convergence of water within the soil column, which can generate stormflow, 
is analogous to infiltration-excess (Hortonian) runoff production at the surface. These 
mechanisms are similar in that water accumulates from above faster than it can permeate 
a soil layer below, and instead must be discharged laterally. Both processes are 
represented in our model — full saturation of the top layer not only generates stormflow in 
layer one but also allows for Hortonian overland flow of any continued excess 
precipitation. Our implementation of stormflow thus smoothly merges with the existing 
mechanism for generation of infiltration-excess runoff (see Stieglitz et al., 1997) and 
should be applicable to areas dominated by Hortonian runoff, such as the U.S. Southwest. 

Evaluation of the Modified Soil Column Profile (MSCP) 

Application of the decrease in both porosity and field capacity with depth in the 
soil column produces a more responsive water table and thus intensifies groundwater 
flow discharge. The MSCP modifications impact different aspects of the model 
simulation than stormflow. Response at the quickflow time scale is not engendered; 
rather, baseflow becomes more reactive. Instead of a redistribution of the timing of 
discharge, the magnitude of response increases. Large storm events can be better 
simulated due to a more responsive water table, and transition recessions — from wet to 
dry conditions — are improved by this model modification. However, unlike simulation 
with stormflow, few additional storm events are captured with the MSCP, and early 
spring wetting from snowmelt is not depicted. In addition, during drier conditions, the 
model still underrepresents discharge between storm events. 
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Our sensitivity analyses show that differences in catchment response to the MSCP 
exist. The steeper catchment of Black Rock with its shallower upslope soils requires a 
greater rate of decrease in porosity and field capacity to best simulate discharge. This 
difference in mean profile is warranted in the single column TOPMODEL framework. 
Because the MSCP is a one-dimensional parameterization of structural variability in 
upland and lowland regions, physical differences among catchments will require different 
parameterizations. The MSCP parameterization does not violate TOPMODEL 
assumptions and thus may be appropriately calibrated. 

Better hydrologic simulation might be achieved with a more precise determination 
of the soil column profiles for porosity and field capacity. The geometric decay functions 
adopted in this study are admittedly somewhat arbitrary. The form of these functions was 
chosen to depict a combination of effects — differences in depth to bedrock among upland 
and lowland areas, as well as bioturbation, and soil compaction— in the absence of 
detailed experimental evidence. In the future, this profile could be more precisely 
matched on a catchment-by-catchment basis to analyses of soil samples. These matched 
functions would need not be monotonic. 

Evaluation of the Joint Stormflow/MSCP Application 

Combined implementation of both stormflow and the MSCP produces the greatest 
improvement in model simulation of catchment discharge. In these simulations, model 
response to storm events during dry conditions better matches measured discharge 
intensities and durations. Water that would have otherwise been lost to the deep water 
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table is discharged via stormflow, and the MSCP intensifies this response for larger 
events. The magnitude and timing of spring snowmelt response is also much improved. 
Stormflow responds quickly to these wetting events, and MSCP permits faster 
groundwater flow reaction to persistent wetting such as during snowmelt. These two 
effects amplify one another. Finally, recession curves during the transition from wet to 
dry conditions are also better represented. The more dynamic water table, produced by 
the MSCP, effects a more realistic groundwater response during these periods. 

Sensitivity analyses also show that different optimum MSCPs exist for the two 
catchments. As for the model runs with the MSCP alone, this variability reflects 
differences in catchment structure and relative soil depths among upland and lowland 
areas. 

Simulations with altered soil column layer resolutions demonstrate that the 
impacts of stormflow and the MSCP are robust within the constraints of the model 
layering scheme. We do not suggest that the stormflow and MSCP methodologies are 
entirely resolution independent; however, provided there is sufficient discretization of the 
soil column, stormflow and groundwater flow discharge simulation appear to be 
relatively insensitive to changes in layer thicknesses. 

Summary 

This study has presented two new modeling strategies: a methodology for 
representing subsurface stormflow; and a modified soil column profile (MSCP) that 
produces a more responsive water table. Both stormflow and the MSCP are adopted in a 
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manner consistent with TOPMODEL assumptions. The stormflow methodology allows 
for discharge from regions of near saturation that develop in the vadose zone during 
storm events. Determination of the stormflow contribution to discharge is made using the 
equation for groundwater flow, and is a function of catchment topography and local 
hydraulic conductivity at the depth of such regions of near saturation. Stormflow 
simulation produces discharge at quickflow time scales. 

The MSCP represents the different soil column profiles of porosity and field 
capacity in upland and lowland regions of a watershed in the single column TOPMODEL 
framework. This parameterization of physical structure creates a more dynamic water 
table and more responsive groundwater flow. In applying the MSCP, large storm events 
are better captured and transition recessions from wet to dry conditions are better 
simulated. Steeper catchments with shallower upland soils appear to require a greater 
reduction of porosity and field capacity to better simulate large storm event discharge. 

The new modeling strategies have been applied to two experimental watershed in 
the northeastern United States. Used jointly, stormflow and the MSCP provide a more 
accurate representation of the time scales of catchment response to storm events and a 
more complete depiction of hydrologic activity. Storm events large and small are better 
simulated, and some of the biases previously evident in TOPMODEL simulations are 
reduced. 
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Figure Captions 

Figure 1. Schematic of the hydrology model. The model couples the analytic form of 
TOPMODEL equations within a discretized column framework. From an update of the 
mean water table depth, TOPMODEL equations and DEM data are used to generate 
groundwater flow and the saturated fraction of the watershed. 

Figure 2. Schematic of the stormflow and MSCP modifications, a) Baseline model soil 
column with a fixed level of porosity and field capacity with depth. Porosity and field 
capacity are represented jointly by the open areas of the column; solid soil is represented 
by the hatched brick filling the left of the column; in this version vertical infiltration 
recharges the water table raising groundwater flow, b) Model with stormflow — gravity 
directs a component of the recharge laterally generating stormflow; c) Model with 
stormflow and the MSCP — as for b) but with reduced porosity and field capacity with 
depth (represented by the expanding hatched brick). 

Figure 3. Calculation of stormflow. a) Calculation of stormflow in layer i begins if 
saturation exceeds 70% of field capacity, b) Water in excess of 70% saturation is 
redistributed in a fully saturated band from the base of layer i upward; all layers below 
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1973 Sleepers River Catchment Simulations 
Figure 5a 

Baseline Model Hvdrograph 



Figure 5b 
Model Hydrograph with Stormflow 
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1999 Black Rock Simulations 
Figure 6a 

Baseline Model Hydrograph 



Figure 6b 

Model Hydrograph with Stormflow 
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Figure 10a 

Spectral Analysis of Sleepers River Hvdrographs 
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Figure 10b - Baseline Model Simulation 
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Figure 10c -Model Simulation with Stormflow 
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Figure lOd - Model Simulation with the MSCP 
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